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ABSTRACT 


Three probability models for defense against mass missile attacks are 
developed. Each model corresponds to a different level of information. Different 
strategies are made possible by the information available. The first strategy, Invisible 
Kill, assigns equal times to the destruction of each missile that is selected for 
attention. The second, a Threshold strategy, assigns a maximum threshold time to 
the same task. The last, based on Perfect Task Information, engages missiles in 
ascending order of their destruction times; the order being assumed known in 
advance. The Expectation of the number of missiles killed is the Measure of 
Effectiveness used to evaluate the models. Numerical results are analyzed through 
simulation. Different strategies are compared showing the effects of information in 


defense effectiveness. 


THESIS DISCLAIMER 


The reader is cautioned that computer programs developed in this research may 
not have been exercised for all cases of interest. While every effort has been made, 
within the time available, to ensure that the programs are free of computational and 
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without additional verification 1s at the risk of the user. 
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I. INTRODUCTION 


Information is known to be an important factor in the determination of the 
outcome of military combat. Consequently, the addition of information-gathering and 
interpretation assets must be viewed as in legitimate competition with the acquisition 
of weapons. In fact, modern "smart" weapons possess within themselves a 
combination of brawn: destructive power (e.g. high explosive warhead); and brain: 
navigation, identification, guidance and homing power (e.g. on-board sensors and 
logic). External C? systems are needed to guide overall decision-making: on the 
defender's side to detect attacks and to schedule defense response; on the attacker's 
side to first soften or confuse (EW measures, deception) the defender prior to a 
focused destructive action against sensitive assets. 

This thesis seeks to investigate the effects and value of information availability 
in a defense against a single mass (simultaneous) missile attack. The analysis is 
carried out from the point of view of a defender operating within constraints of a 
fixed time window, presumably until destruction of the attacking missiles is no 
longer possible at the level considered; "leaking" missiles may be dealt with at an 
ensuing defense level. For simplicity we postulate that the objective of the defender 
is to maximize the expected number of missiles destroyed at the level under 
consideration. The defender is using some weapon system (interceptor missiles, 
irradiation. weapons,..). This generic weapon system is characterized by its 
probability of acquiring a target, acquisition or setup time, and rate of kill. These 
parameters are supposedly design features intrinsic to the system. It is possible that 


trade-offs between the weapon system characteristics (parameters) exist, in which 


case our models may provide some insight in deciding the desirable parameter 
configuration. 

The three models we develop are of a probabilistic nature, each of them 
corresponding to a different information level available to the defender. A higher 
information level should enable the user to enhance his strategy and consequently 
achieve a gain in defense effectiveness. The general methodology followed is to start 
with the basic assumptions of the model and strategy, and work our way analytically 
to explicit expressions for the chosen measures of effectiveness. Numerical results 
of the probability models are checked against those of discrete-event simulations. 

It must be noted that, in this context, the term "information" is used in a 
specific but informal sense appropriate to the present application. There is no direct 
connection with such classical concepts as Shannon information or Fisher 
information. Finally, note that the present problem may appear in other guises, in 
which a number of tasks must be performed under time constraints; e.g., in repair 


and maintenance or flexible manufacturing. 


Il. MASS ATTACK WITH INVISIBLE KILL 


A. INTRODUCTION 

This chapter deals with the following situation: A defender is under a mass or 
simultaneous attack from a possibly large number of incoming missiles (a cluster’). 
Available to the defender is an unspecified outer-layer weapon system, D1, which 
enables him to shoot at least some of the missiles. Each missile, upon being detected 
must be acquired, and only then can one attempt to shoot it down. The weapon 
system only allows the defender to engage one target at a time, so he faces the 
problem of how to optimally assign time to each missile so as to maximize some 
measure of the number of missiles destroyed. It is understood that one is working 
in a restricted time window: if not destroyed, the missiles will proceed towards their 
targets, passing through a secondary defence level. The missiles that are not killed 
by Di are said to leak to the next layer. The general objective is to minimize 


leakage. 


B. ASSUMPTIONS 
In the probability model the following simplifying assumptions will be used: 
e The number of attacking missiles is known. 


* The probability of acquiring each target 1s known and fixed; range dependence 
is ignored. 


* Given that the missile is acquired, the probability distribution of the time to 
destroy it 1s known. 


* It is not known to the defender whether a given missile was successfully 
acquired or, if acquired was actually killed. This is the minimal information 
situation; we call it Invisible Kill. 


e There is a fixed setup time expended each time a new attempt is made to 
acquire another missile. 
Finally, the strategy of the defender will be to allocate time evenly between those 
missiles selected for attention by D1. Not all attacking missiles will be chosen under 
optimal conditions; some will be allowed to leak, presumably with warning 


transmitted to the second (inner) defence layer. 


C. MODEL 

The following notation will be used: Let n be the number of incoming missiles. 
The probability of acquiring a target is p , and 4 - 1-p the probability of 
acquisition failure. Denote by X the random time to destroy a missile which has 
been acquired. T is the time available for destroying the targets - the duration of 
the window of opportunity. There is a fixed time necessary for each acquisition 
which will be denoted by w. 


It is assumed that the distribution of X is known. Call its CDF F,(x), where: 


X ~ F,(x) = Prob{ X < x | Target was acquired } . (2.1) 


Suppose one chooses to engage a number, m, of the incoming 
missiles; m < n. Recalling the declared strategy of assigning equal time to each 
one, then Т = mw+mx; here x is the time allowed to destroy a single missile. 


It follows that: 


- уу. (22) 


Since each target must be destroyed in less than time x, Prob{destroying a 
target} = Prob{a target is acquired upon attempt} * Prob{a target is destroyed in the 


time available: max(0,7/m-w)  }; Le.: 
| ] Т , 
Prob (destroying a target) = p . F.( — -w). (20) 
т 


Note that since X > O with probability one, automatically Ё, (х) = 0 if x < O. 
If it is assumed that each attempt to acquire and destroy a target is an 
independent Bernoulli trial, then the number of targets destroyed out 


of m, N(T), will be distributed 
| Г 5 
N(T) ~ Binomial(m ,pF,(—-w)), (2) 
m 
and the expected number of targets destroyed will be 
Е Т > 
E[(N(T)] = mpF,(—-w). (2.5) 
т 


The current problem is to find the optimal value of m, i.e., one that maximizes this 
expected value. Other criteria and controls are of course possible. Note that if the 
targets are known to be different, then a strategy that schedules different times to 
each type would be reasonable. An extension of this formulation might include a 
preliminary classification of target type and subsequent allocation of effort to the 


most profitable types first. 


D. MODEL DEVELOPMENT 
The optimal value of m will now be derived: namely the number of missiles 
the defender should engage in order to destroy the maximal expected number of 
them at his particular defense layer. 
The expected number of destroyed or killed missiles is given by (2.5). We 
wish to find the maximum of this function. Set the derivative SENN 0 


дт 
implying that if the optimal value, m = m°, is interior, O<m°<n, then 
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Flz-w) =2f(z-w), (2.6) 


where the parameter z=T/mº. In principle there can be several such stationary 
points, and furthermore the optimal number can be a boundary value: m° =n. 
In particular suppose X ~ Exponential(à), where the rate parameter may be 


interpreted as the number of acquired missiles the weapon kills per unit time. Then: 
F, (x) = 1 -e** x>0; (2.7) 
The expected number of destroyed missiles becomes: 


E[N(T)] = трП - е Жш. (2.8) 


To maximize this function of m, take the derivative and set it equal to zero. 


After some algebra, get: 


In(m) - In(m+aT) + 24 - Aw = 0; (2.9) 
m 
An approximate explicit solution for m”, a stationary point, even if somewhat 
inaccurate, would be desirable to show its dependency on the various parameters. 


Start by rewriting (2.9), and using a Taylor expansion: 


о е 
т 
(2.10) 
2 
РЕГ ae 
т т 2m? 





solving for m, yields the approximate explicit result: 


m =T ES (ТІ) 
2w 


This number will be a positive real. We want an integer solution; therefore the 
greatest integer not exceeding m’, [m'], and the next lowest integer, <m >, 
must be checked to see which one has the largest objective function (2.8). Call this 


* 


number m`’, then the approximate optimal m? will be: 


0 


m =min(n,m'). (2.12) 


Accordingly the corresponding optimal value of the expected number of 


missiles killed is, again, approximately: 


т 


© ы, 2.1 
E°(N(T)] = т°р [1-е Б |, бЭ 


Note that the second derivative of E[N(T)] 15 
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SEND] . -p (АТ) — (2.14) 
дт? m^ 


therefore 


9* E[N(T )] 


5 <0, Ут. (2915) 


дт 


The function E[N(T)] is concave throughout the domain of m. In this 
case m°, if obtained as described above, will be the approximate unique global 


maximizer. 


E. SIMULATION 

To verify the analytical results obtained thus far, a terminating simulation 
[Ref. 1] was performed. The measure of effectiveness estimated is the Expected 
Number of Missiles Killed per Attack. The attacks were replicated independently, 
15000 times, with several combinations of input parameters. The FORTRAN codes 
pertaining to the Event routines are displayed in Appendix A. For Timing and 
Calendar manipulation routines the SIMUTIL package [Ref. 2] was employed, as 
was the LLRANDOMII Random Number Generator [Ref. 3]. 

The results of the simulation were written to CMS files, and used to generate 
several plots of response surfaces. One of these plots is shown in Figure 2.1, with 
E[Missiles Killed] as the response variable. 

For example, if one takes the time window to be 700, 20 attacking missiles, 
the rate of kill 0.01, the acquisition time 10 and the probability of acquisition 0.85; 
then solving (2.11) yields m* = 15.652476. Checking m* = 15, 16 with (2.13) 
reveals that m "' 2 15; therefore m? - min( 15,20) = 15. It can be observed from 
Figure 2.1, with these values, that the approximate analytical optimal seems to be 
in good accordance with the simulation results. Note that obtaining the precisely 
correct value of m is not important; the value of the expected number killed is not 


very sensitive to that value. 


INVISIBLE KILL 
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Figure 2.1 Simulation Results 


Е. USE OF DECOYS BY ATTACKER 

It might be of some practical interest to discuss the influence of the use of 
decoys, in the results of the previous strategy for Invisible Kill. By assumption the 
decoys will be randomly mixed with the missiles within the cluster. As information 
available is minimal, the decoys will be undistinguishable from the missiles for the 
defender; as a result the defender employs the same strategy regardless of their 
presence. 


Let r be the fraction of targets which are true missiles, Le., 


= number of missiles (2.16) 
number of missiles + number of decoys 
Using the current independence assumptions, the number of missiles destroyed 


in m trials will be distributed Binomial and its mean will be: 
и Т 
E [missiles destroyed] = mprF,(— - w). (2980) 
т 


If the destruction time is again taken to be Exponential, the same steps already 
taken show that the optimal number of targets to engage does not vary with r. The 
expected number of missiles destroyed using the optimal strategy decreases, being 
multiplied by a factor of r. This result is rather intuitive. Since there is no way the 
defender can assess the composition of the attacking target cluster, defense 
effectiveness depends directly on the fraction of decoys. Eventually an attack with 
large numbers of decoys will succeed in diluting the defender's effort, thus 
achieving a high rate of penetration to a subsequent defense layer. 

It is interesting to conjecture the payoff from a subsystem capable of 
classifying attacking entities as either missiles or decoys with some probability 
("skill"). It may be that tradeoff between some of the other parameters, 
e.g., T,w, or A, and classification skill can be beneficial for the overall system 
effectiveness. This tradeoff is susceptible to quantitative study conducted in the 


mode described throughout this thesis. 
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Ill. VISIBLE KILL - A THRESHOLD STRATEGY 


A. INTRODUCTION 

The former chapter dealt with a situation in which there is no direct 
information pertaining to when, or if, a given target is destroyed. The previous 
assumptions are next modified in a way such that the defender recognizes the instant 
when a target is destroyed (VISIBLE KILL), if such an event occurs. Upon 
destruction of a target the defender immediately attempts to acquire another (if any 
are left). It should be noted that acquisition is still ‘not visible’, i.e. after the fixed 
acquisition time w, the defender does not know if the acquisition was a success. 
Thus D1 may have to wait indefinitely for news that a kill has occurred, which leads 


to consideration of a time threshold strategy or policy. 


В. THRESHOLD MODEL 

A threshold policy is adopted: the defender will allow a maximum threshold 
time ‚т, to be allocated to destruction of a target. If, after time т elapses, the 
target 1s not observed to be destroyed the defender will go back to acquisition mode. 
Such a policy is easy to implement and can be justified by dynamic programming 
under certain circumstances. 

Define Z(:), the total random time to kill one target with a threshold equal 


to t. Then, from the assumptions made, the following recursion must hold: 


wt+x xx 
Z(t) |X=x = сш 


wo + (т Ж к 


where X is the random time to destroy an acquired target, in a continuous 
engagement. Z(t) is a random variable having the same distribution as Z(t); it 
is the time to shoot down a target if the first attempt fails and the process of 
acquisition must start over. The CDF of X is known to be Ё, (х). 

The upper branch of (3.1) occurs when acquisition is successful and the target 
is shot within the prescribed threshold. This has probability pdF,(x). The lower 
branch occurs if either acquisition fails, or acquisition succeeds but the missile is not 
shot in the allowed threshold time. This has probability q +pF,(t). Thus we 


obtain the unconditional expectation: 


w+ [pxdF,(x) +т[д +рЕү(т)] (3.2) 
2 


p Fy(t) 
where F (x) = 1 - F,(x). Integration by parts in the numerator gives this slightly 


simpler form: 


w +p | F dx + та (3.3) 


E[Z(x)] = "d 
X 


The defender wishes to maximize E[N(T)], with x being the decision 
variable. Therefore we must express E[N(T)] in terms of the expectation (3.3) just 


derived. 
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If we interpret each missile destruction as a point event, then Z(z) would be 
the random interarrival time. If we had a ‘continuous supply" of missiles, and a 


large ume window ,T, then from the elementary renewal theorem [Ref. 4]: 





13 
EINT) =- (3.4) 
Е(Жт)) 
In fact, a correction term might be added [Ref. 4): 
_p2 
Е[М(Т)) ~ T ‚ Гав) ЕЕ [200]. (3.5) 


E[Z(x)] 2E [Z] 
This requires the computation of the variance of Z(t) which is more difficult and 


leads to a complex expression. In this paper we retain the approximation (3.4), 


accepting the consequent loss of accuracy. Thus we use 


T p F,(x) 
ЕМТ)] = - 
а (3.6) 
w +p | FyG)dx +79 
0 
Тһе defender's problem is to find: 
т’ = агетах Е МТ)] (3.7) 


De qa T 
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т°, the optimal threshold can be looked for among the critical points on the 


interior of its domain, or at t =T. 


For example, say X ~ Exponential(à). After substitution of the appropriate 


expressions for the distribution: 
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о ЗЛЕ 
ВАТ = =”; 


(3.8) 
wri (1 -e ^')««q 


Simulauons were conducted to study the dependence of E[N(T)] on the 
parameters À, w, t, p. Plots resulting from these simulations are displayed and 
discussed on Appendix B. 

Again we seek the critical points of E[N(T)] by finding the zeros of its 


derivative: 


DEINO ОЕ ЕА (3.9) 
OT q 


or equivalently, by taking logarithms of both sides: 


BRE CE (3.10) 
q 


If p - 1, it can be verified that x? - T, otherwise the optimal threshold is to be 
found as a solution to (3.9) or (3.10). 
Two approximate solutions for (3.9) are provided. The first is obtained by 


rewriting (3.9) as 
Ат Ам 
e^ sme 8 (3.11) 
and expanding the exponential term: 


ER DU (3.12) 
2 


Substituting (3.12) into (3.11) and solving for t yields: 
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25-225 (3.13) 


т = — 


Aq 
A second, more accurate solution, is derived by taking (3.13) as the initial 


guess and performing one iteration of Newton’s method to solve (3.11), resulting in: 


2wa 


| Ba a) ШО 
o q q 


ЖО e 


Optimal threshold values obtained with the approximations (3.13) and (3.14) 


(3.14) 
T 


are shown in Table 3.1, together with numerical solution of (3.10). 


i — —- 3.1] OPTIMAL THRESHOLD APPROXIMATIONS (pz0.85) 


E vs | 
18.164104 | 
—— — — 
— 14) 20.659227 
G. —— 
в اف‎ 
29.264013 
Г ем." 10) 128.640050 | 45.983242 28.435654 
20.0 51.639778 
42.703706 


Although the approximations are not always accurate, they are still useful since 
the expected number of missiles killed is not very sensitive to the threshold value 


(see Appendix B). Table 3.2 shows values of E[N(T)], obtained from (3.8) with 


optimal thresholds given by (3.10), (3.13) and (3.14). As it is seen, the three 


methods of deriving the optimal threshold produced very simular results. 


Table 3.2 EXPECTATION OF MISSILES KILLED (T=700, p=0.85) 


— ом | оз | oa _ 


20.0 


" — 
— 14) 
810... 10) 
"WU u 
(3. в 
(3. | сю 
(3.13) 
(3.14) 
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С. SIMULATION 

Optimal values, E?[N(T)], obtained from (3.8) and (3.10), were compared 
to those of simulations, for several values of A. The simulation code is very 
similar to the program "Threshold Simulation" given in Appendix A. The results are 


displayed in Figure 3.1. 


Visible Kill 
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Figure 3.1 Threshold Policy - Results 





As can be seen, the model slightly underestimates the simulated values, 


although its adequacy is seen to improve for larger values of the killing rate. 
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IV. DEFENSE WITH PERFECT TASK INFORMATION 


A. PROBLEM STATEMENT 

Suppose one is given n tasks, with the objective of completing the maximum 
possible number of them within a fixed time window of length T. (E.g., a batch of 
jobs arriving at a repair facility at the beginning of the day; a cluster of missiles that 
must be destroyed in limited time). Assume, as an example, that the duration of each 
task is Known, and drawn from an Exponential distribution. Since the objective is 
to maximize the number of tasks completed, one obvious strategy is to process the 
tasks in ascending order of their completion times. This policy will be referred to 
as a Smallest-First (SF) rule, in the same fashion Coffman [Ref. 5] uses Largest- 
First. In the context of our missile defense problem the above formulation is 
artificial because task times would never be known in advance. Note, though, that 
application of the rule only requires knowledge of the order, and not the actual task 
(shoot-down) time. Therefore the results provide an upper bound on the 
effectiveness of the outer defense layer, D1. 

Denote by X, the time for completion of the i" task. These times are 
independent and identically distributed exponentially with parameter À, their i“ 


order statistic being X. If the SF rule is employed, then for some j < n: 


D' 


X 


^ (4.1) 


+... + Хр 5 * Ау, 


in which case only j tasks may be completed or "packed" in time 7. Call the 
number of tasks packed with a SF rule N(T), a random variable. The aim of this 


chapter is to find the probability mass function (PMF) of N(T), or accurate 
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approximations to that PMF. This can then be used to evaluate measures of system 
effectiveness, which can be compared to the results of applying other rules that are 
based on less information. 

Similar problems have been treated in recent literature, although in a more 
general way and with different objectives and methodology than that presented here. 
For example, Coffman, Fayolle, Jacquet and Robert [Ref. 5] obtained asymptotic 
results for the expected number packed, under a Largest-First (LF) rule; with the 
samples being drawn from a U(0,1) distribution. Coffman, Flatto and Weber 
[Ref. 6] derived the asymptotic value of the expected number of intervals selected 
with an optimal ‘threshold’ rule, for a class of distributions. Theirs is an "ОМ- 
LINE" rule, i.e., according to their own definition, the decision to pack or not the 
interval (job) is made at the time of inspection. The missile problem, 1.e. of rank- 


ordering task times, would also necessarily have to be solved "ON-LINE". 


В. EXACT PMF OF МТ) 


The following well-known results are relevant to the derivation of the PMF 
of N(T). 
1. Distribution of Spacings 
The distribution of the differences between exponential order statistics 1s 
next presented. More detail can be found, e.g.. in Barlow and Proschan [Ref. 7]. 
Let A ул DE 11.4. ехропеппа! with parameter A. Define the 


spacings, the differences between successive order Statistics, as: 
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The normalized spacings nD,,(n-1)D,,..,D, are again i.i.d. exponential with 
rate à . Let these normalized spacings be denoted by Y,,..,¥_, respectively. 


They can be expressed in terms of the order statistics in the following way: 


= ПХ р 


Y, = (п-1)Х„, - (п-1)Х\, (4.3) 
ТР (n -]* 1) Xj = (n -J +1) Xj- Р 


The last equations, if solved for the order statistics, yield: 


_ 1 
Ха E _1, 





Xo = +Y, + >Y, (4.4) 





1 1 1 
Х ES Е + ب‎ С a ecc T Т • 
9) n | n-1 ? п-ј+1 1 


Define the summation of the order statistics: 
- _ QS j-itl 
т № Ха - 23 Y 


i? 


iz] iz] n-1+1 (4.5) 





Taking the expectation and variance of the right-hand side of (4.5) produces: 
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ГЕ Шеге 
Eis) === 
i ‚ п-із1 
(4.6) 
г 
Var(S] = — 
E з 


2. Convolution of Exponentials 


The following result, for the convolution of independent Exponential 
densities, will be used in deriving the distribution of 5. ЕЕ 1 


ЕТЕ X, have densities 


fy, Ee TOI (4.7) 


Then the sum À, +... + X, has density: 


J J 
f = П 3» в,« = 
: t=] i=l | 


with the B, ; defined by the following product: 


J 
Bue (4.9) 
К=] 


kei 
3. Derivation of the PMF of N(T) 


In the current case Si as defined in (4.5), is a linear combination of 


j i.i.d. exponential distributions, as seen in section [ILB.1. Also by a simple 


transformation of scale: 
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усу ~ Exponential 7%, (4.10) 
n-i Jir 


thus 5 can be viewed as a sum of j independent exponential random variables 


with different rates. Then applying result (4.8) yields the distribution, after some 





algebra: 
J ا‎ 
FO - Н)! Ж. خا ا‎ тім е", 
Jli IJK j-t+l (4.11) 
kzi 
ана A EUR MO 
To derive the PMF of N(T) the density above is used: 
оО =J) = Еу (D - Fs) ; (4.12) 
with F(T) = P(S > T0) -IL-EFSIT), 
/ ) ) 


where Е; (0 is the CDF of S, and Е, (0 is obtained from the density (4.11) by 


integration: 


Е,(Т) = J Һа: 


- (95 (f Де пу aci ііі АТ. (4.13) 
= Ке 7-1 — 


defined forj =1,..,m-1. 


Using (4.12) we have the desired PMF for 7 =1,...,n-2. The 


values 7 =0,n-1,n have not been covered so far. 
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POG) any NEUES (4.14) 
i=] í 


since S, is a sum of i.i.d. exponential Random Variables, then: 


S, баппат л) A4 >Q; 
(4-15) 
= РМТ) = 





>. 
=] 
aM: 
| 
حر‎ 


Em 
If no task is completed, the shortest time to complete a task must be 


bigger that the time window, 1.е.: 


РАМТ) =0) = Е; (Т) = е "АТ. (4.16) 


Then by a generalization of the domain of definition in (4.13), 


0 f j <0 
(4.13) ЕІ-і.,..п-1 
r S (4.17) 
FS) ET Tan “(AT — 
jo j! 
1 if j2n-*l 


With this broader definition of F,(T) (4.12) holds for the whole range of 
J 
j,{0,1,..,n}, and an exact PMF of N(T) has been derived. 


C. APPROXIMATION TO THE PMF OF N(T) 
To recapitulate, an exact but complicated formula for the PMF of N(T) has 


been derived earlier; see (4.17). We now desire an approximation, which might be 
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easier to use. This can be achieved via the distribution (CDF) of S j in the following 


way: 


PINTY =J) = RNO al RD) 


(4.18) 
- Е,(Т) - Б, (Т). 


1. Gamma Approximation 

Analysis of simulation data, e.g. using probability plots, suggest that 
the 5, are approximately distributed Gamma (see Appendix C for details). If the 
approximation is good enough it can be used to compute values 
of Р(М(Т)-/). The values of the Expectation and Variance of the 5, are 
available in (4.6). 

If 5, is to be distributed Gamma with parameters shape а, апа 
rate В, then: 


Q. 
EIS) = апі Var[S] - — ; (4.19) 


J 


ое 


from which the parameter values of the approximating distribution can be derived, 


since the values of the expectation and variance are known: 
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(4.20) 








It was seen that this approximation produces good results in computing the PMF of MT) 
using (4.18). Another approximation method is presented next. 
2. Normal Approximation 
If we assume that 5, can be represented by a Normal distribution then 


the resulting Normal approximation,with ®() denoting the probability mass in the 


left tail of a Standard Normal, is: 


ОС =e 


(4.21) 





Т-Е[5 ] 4 I ES i | 


[Varis] /VarlS, 1] 


for j-1,..,n-1; and with EIS, Varls |] defined by (4.6). 








D. RESULTS 

Once more analytical results were checked through a Monte Carlo simulation. 
The simulation program replicated independent “packing experiments". The values 
of "Number Packed" (or tasks completed) thus obtained, were then used in an APL 
function to compute the “experimental” distribution of N(T). This distribution is 
plotted together with the corresponding approximations and exact PMF, (4.12). The 
inclusion of a continuity correction in the Normal Approximation provided no 
substantial improvement. The Normal approximation was seen to work better than 


portrayed in Figure 4.1, for cases with bigger values of n, as should be expected. 


The Gamma approximation produced results quite similar to those of the exact 


formula. 
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Figure 4.1 Simulation and Analytical Results 


Suppose an unsophisticated policy is used, so jobs are taken randomly, without 
ordering. In that case МТ) ~ Potsson(A) , and E[N(T)] = AT. The PMF of 
the number of jobs completed in this way, was compared with the corresponding 
PMF under a SF rule. This gives a measure of the benefits of employing the latter 
policy. Results are summarized in Table 4.1. It can be seen in the bottom row that, 
for this particular set of parameters, the SF rule accounts for an increase in Expected 


value of about 85%. 
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Table 4.1 COMPARISON OF RESULTS 
FOR THE — — OF TITHE NUMBER OF TASKS COMPLETED 


LERRA ont SIMULAT | NORMAL | 
| - [9 [0030197 
— — 0.105691 
Го | о 10000031 | 0.000056 | 0.184959 
[4 | 0.037800 | 0019231 | 0035266 | 0.035456 | 0.188812 - 
| € | 0.152300 | 0195998 | 0.158881 | 0.156441 | 0.132169 | 
[6 | 0315000 | 0.342860 | 0.310115 | 0.309978 | 0.077098 | 
И 0.152500 _ 0.122254 


0.035700 0.047305 0.036936 0.036370 0.006559 | 


0.002600 fe 0.003315 а 
| avG || 6.487500 | 6.284765 | 6.458033 | 6.493038 | 3.500000 | 
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Е. APPLICATION TO MISSILE DEFENSE 
We now proceed to incorporate the results obtained from the packing problem 
to a missile attack situation. As before, the defender must acquire and shoot the 
maximum number possible of n incoming missiles, within a restricted time window. 
1. Assumptions 
The missiles are acquired with probability p, acquisitions being 
independent Bernoulli trials. A fixed time, w, is spent for each acquisition. The 
defender, knowing the ordering of the "missile killing times”, will select them for 
acquisition with a SF rule. It is not recognized if acquisitions are successful or not. 
Consequently, the first time the defender fails to acquire a missile the residual time 


window is expended and no more missiles are shot. Presumably the performance 


2 


could be improved by introducing a threshold, as was done earlier, but this option 
has not been studied. 
2. Model 
Let N(T) now represent the number of missiles killed in time 7, under 


the above assumptions; then, 
P(M(T) 2 j) = p/F,(T - jw); (4.22) 
resulting in the PMF: 
P(N) =) = pF E - p"! F, (T-G*Dw). )4.23( 


Recalling that F, (1) = 1 -Е; (0) сап be obtained from (4.13). 

The result (4.23) applied to a given set of parameters is plotted in Figure 
4.2 together with corresponding values obtained through simulation. 

For the purpose of comparison the results of the model were checked 
against those obtained from employing the unsophisticated strategy, i.e., taking the 
missiles in random order. In this case the distribution of N(T) is obtained in the 
same fashion, by using (4.23), but with the difference that now 5, is distributed 
Gamma with shape j апа rate A. The model is seen to agree reasonably well 
with the "experimental simulation data" as further detailed in Table 4.2. It is also 
verified a substantial improvement by using the SF strategy. 

A noticeable feature of the distribution is that it is no longer unimodal, unlike 
similar distributions for the packing problem. In plots generated with different rate 
values one local minimum appeared between two maxima, one of them located at 
zero. Figure 4.3 depicts a situation in which the rate is 0.1 and 20 missiles are 


present, all other parameters remaining equal to the previous situation. As can be 
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seen the minimum is now more sharply defined, and the absolute maximum is at 


Zero. 


Perfect Task Information 
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Figure 4.2 Model and Simulation Results 


Table 4.2 PERFECT TASK INFORMATION PMF 


SIMULATION | MODEL NO STRATEGY 


0.148574 0.182953 | 
0.117059 0.238861 
0.113057 0.259735 


0.205103 0.028368 

0.107554 0.005541 

0.000500 
0.000011 0.000002 
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V. CONCLUSIONS 


A. MODEL VERIFICATION 


In the preceding chapters we developed three models for defense against 


missile attacks. For each model analytic results were compared with simulations, for 


some range or combination of parameters. To summarize, our three models 


correspond to the following situations, in ascending order of information availability: 


Invisible Kill 


Unknown: If acquisition successful. 

If and when missiles shot. 
Strategy: Allocate time evenly between all engaged missiles. 
Decision: How many missiles to engage. 


Visible Kill 

Unknown: If acquisition successful. 

Known: Time when missile killed (if killed). 

Strategy: Allocate a max threshold time to each missile engagement. 
Decision: Length of threshold. 

Perfect Task Information 

Unknown: If acquisition successful. 

Known: Times of destruction of each missile. 

Strategy: Engage missiles with SF rule. 


Throughout our work, whenever numerical results were necessary, assumptions 


on the values of the parameters had to be taken. We have used what seemed to be 


reasonable" values of killing rates, acquisition times, time window, etc. These, 


however, may differ from actual parameter values to be found in existing weapon 


systems. For this reason we now seek to "validate" the models in more general 


conditions. We wish to verify and compare our models in a reasonable number of . 


parameter (factors) combination without too much simulation effort. The response 
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vanable to consider for validation/comparison purposes will be the expected number 
of missiles destroyed. The models depend on five input parameters or 
factors: n,p,w,A,T. The same distribution for the Killing Time will be used - an 
Exponential. In order to generate input parameter values we chose a two level 
fractional factorial design, ا‎ [Ref. 9], which seems to attain the above stated 


objectives. This design corresponds to the input data set displayed on Table 5.1. 


Table 5.1 INPUT DATA 


NEL Hs O. анар oe e ты, 





005 | 0 | 0%. 200.0 


Le Le pae aa 70.0 
700 | 


0.05 10 700 | 


0.05 50 5.0 200.0 
0.1 2000 — | 
m so | 10 700 | 
| | 005 | .s0. M "oss '4 150 2001 


1. Invisible Kill 
The simulation was carried out in the following way: For each run, the 
optimal value of the number of missiles to engage was computed according to the 
model equation (2.11). In each run the simulation replicates enough independent 
attack experiments to get estimates of the average number of missiles killed and 
95% confidence bounds. The same average was computed from the model equation 


(2.13). Results are compared in Table 5.2. 


Table 5.2 INVISIBLE KILL VERIFICATION 


— — 0 1.264240 
7.768698 
2.981409 
TU 
6.658032 
2.723549 


2.660858 | 2.682709 2.704558 2.675591 — 
















г 
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As further depicted in Figure 5.1, all runs produced satisfying results. 


Considering that (2.13) is an approximation, the model values seem to be in 


accordance with those of the simulation. 







E [Missiles Shot] 
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Elqure 5.1 Іпуізіріе Kill Verification 


2. Visible Kill - Threshold Strategy 
We first determined the optimal threshold for each run by numerical 
solution of (3.10), (when p = 1,q =0 the optimal threshold will be T). The 
simulation is executed for each combination of factors (run) and optimal threshold. 
Analytical values (model) are obtained from (3.8). As can be seen in the following 
table and figure, the results obtained from the model are not always accurate as 


compared to the simulation. In fact the model predicts an expectation of about 13 
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missiles destroyed, when only 10 are present in the attacks corresponding to run 2. 
This is an expected result considering the crude approximations taken in Chapter III. 


Otherwise we find the agreement model/simulation to be quite acceptable. 


Table 5.3 VISIBLE KILL VERIFICATION 
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Piqure 5.2 Visible Kill Verification 


3. Perfect Task Information 
Proceeding in an analogous way we compare the expected number of 


missiles shot obtained from simulation, with model results computed with: 


ЕМТ)] = 3, jP(N(T) =j); (5.1) 
із0 


with P(M(T) =j) being given in (4.23). 


The fit between model and simulation is seen to be excellent in this case. 
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Table 5.4 PERFECT TASK INFORMATION VERIFICATION 


LOWER | AVERAGE MODEL 


3.972124 
9.888283 
11.701859 11.743412 
14 
AS | 10613348 


ба ттс i 
| 16 | 4.589483 | 4.664865 | 4.740247 | 4.691323 | 



















B. COMPARISON OF STRATEGIES 

Given that not all the models produced equally reliable results, we chose to 
compare the different strategies using simulation. Again, the input data set of Table 
5.1 is used to generate the values depicted in Figure 5.4. 

Clearly the SF strategy in the Perfect Task Information situation shows the best 
overall performance. This is nota surprising result, since this strategy benefits from | 
the highest level of information. Also the Threshold strategy comes next in overall 


performance, and Invisible Kill produces the worst results. More unexpected is the 
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Figure 5.3 Perfect Task Information Verification 





fact that in three cases (runs 4, 10, 14) Threshold wins from Perfect Task Information, 
which is even outperformed by the modest Invisible Kill in run 14. Although these 
exceptions do not alter the conclusion that models with access to more information 
perform better overall, it is somewhat surprising that this situation does not occur for 
every combination of parameters. A closer look at the conspicuous runs reveals that 
in all of them p =0.85, T - 200. Thus the explanation for the underperformance of 
Perfect Task Information must lie in the fact that acquisition is still invisible (see 
IV.E.1). This implies that, given sufficient time, the defender eventually fails to 
acquire. Unlike with Threshold and Invisible Kill there is no chance of reacquisition. 
However, we see that in more favorable conditions, particularly when acquisition 
occurs with probability one, the SF rule can enhance performance up to five times 


(run 13). 
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Figure 5.4 Comparison of Strategies 


The shortcomings in the Perfect Task Information strategy, can be overcome 
by refining the model assumptions in IV.E.1. The proposed refinement consists on 
including a "threshold" time, such that, after this time has elapsed without target 
destruction being observed, the defender reattempts acquisition of the same target. 
Various thresholds can be used (e.g., the largest task time). The most natural 
thresholds will be the known task times (in this case, knowledge only of the 
ordering is not sufficient). 

With this improvement the Perfect Task Information model clearly dominates 
for all simulated runs, as Figure 5.5 shows. We do not provide analytical results for 


this improved model. 
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Translated into a context where the weapons of the defender consist on 
interceptor missiles, these findings seem to strongly favor a system that enables the 
defender to assign each interceptor according to the shortest time for interception 


(SF). 
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Figure 5.5 Comparison with Improved Perfect Information 
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С. APPLICATION EXAMPLE 

Suppose one is planning a weapon system to face a threat that allows a time 
window T, and that we can forecast that n missiles will be used in each attack. 
We have chosen one of the three models to be the most adequate to describe the 
weapon system that can be built. Then we have, from our analytic work, the 


function: 
E[N(T)] = f(A, p, w; T,n). (5.2) 


The arguments T and m are dictated by the threat, while A,p,w can be 
incorporated in the design of the defense system, within technical feasibility and 
resource constraints, and are thus subject to decision. Further, suppose we can 
estimate the cost functions of increasing the probability of 
acquisition , Cy: increasing the rate of kill ,c,, and decreasing the acquisition 
time ,c,. The decision for the optimal combination of parameter to design into the 


system can be obtained by solving the nonlinear program: 


Мах ИМТ] 
54. c (À) +c (p) +c (w) < R 
ОС Аин (5.3) 
0<р<1 


O< WSW on: 
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D. RECOMMENDATIONS 

The models developed in this thesis may provide a basis for further work 
leading to realistic models of missile attack situations. For this objective to be 
attained, "real world" data and system specific information must be used to validate 
the models. Further investigation may be needed to verify the nature of the 
distribution of the killing times, or the robustness of the models to distributional 
assumptions. 

A possible refinement would be to drop the assumption of invisible acquisition 
which is common to all three models but heavily penalizes the model of Perfect 
Task Information. 

It might also be useful to consider other objective functions. For example, one 
may argue that in a scenario where some missiles are almost certain to leak, it 
makes more sense to attempt to minimize the maximum number of leaking missiles. 
Under some circumstances it may be desirable to maximize the probability of killing 
all attacking missiles; under others it may be useful overall to expend some 
resources to pass information about leaking missiles (course, estimated time of 
arrival) to a secondary defense level, D2. There are numerous other options and 


opportunities. 
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APPENDIX A. FORTRAN CODES 


In this appendix we list the most important programs and routines. The codes 
are written in Fortran77 and were run on an Amdahl 5990-5000, at the NPS 
Computer Center. 

The simulation programs shown were used to verify the numerical results 
obtained with the probability models. The version of the Perfect Task Information 
model presented is the one without the threshold improvement. A program to 
generate exact (analytical) results of the Perfect Task Information is also included. 
Exact computations for the other two models are simple and can be done on a 
programmable calculator. 

All simulations terminate either with the destruction of the last missile or when 


the window elapses. 


A. INVISIBLE KILL 

The main program, "Mass Attack with Invisible Kill", implements the 
simulation to verify the model developed in Chapter II. The program executes 
independent replications of the simulation of a missile attack, to compute the 
expected number of missiles killed. This is done for the 16 sets of parameters 


displayed in Table 5.1, and with т” found as shown in Chapter II. 
Two event routines are present: "Acquisition" and "Shooting. The simulation 
also requires the SIMUTIL package [Ref. 2], LLRANDOMII [Ref. 3], and 


subroutine STAT. 
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A descriptive list of the variables follows: 

"ProbAcq' is the probability of acquisition denoted by p in Chapter П. 
“Time Window” 15 the time available to destroy the targets, Т. 
"IimeAcq' is the fixed acquisition or setup time, w. 

“TimetoKill” is the time allowed to shoot down a single target. 


“Rate” is the killing rate, and also represents the rate parameter of an 
exponential distribution, A. 


"Average", "Variance" and "IntervalWidth" are statistics of the number of 
missiles killed per attack, taken from the sample of all the simulated attacks. 
These statistics are computed by subroutine STAT. 


"Lower" and “Upper” are the 95% confidence bounds for the expectation of 
the number of missiles killed. 


"Iseed" is the seed for the random number generator. It is initialized through 
a "Read" statement at the beginning of the simulation. 


"NrMissiles" is the number of missiles in an attack, n. 


" " 


m" is the optimal number of missiles to engage. It is computed through 
approximation (2.11), as described in Chapter II. It is an input parameter to the 
simulation. 

"Processed" counts missiles as they are generated. 


"NrAcquired" counts the number of successful acquisitions. 


"MissilesDown" counts the number of missiles shot. It is used as an input by 
subroutine "STAT" to compute the average number of missiles killed. 


The integer variables "i" and "Kk" are indices used in "Do" loops. 


"Acquired" is a logical variable and has value .True. when the missile being 
processed has been successfully acquired. 
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Destroyed" is a logical variable. It has value .True. when the missile being 
processed is destroyed. 


"“Tsample" is the random time to destroy a missile sampled from ап 
exponential distribution with rate A. 


"Now" is the current simulation time. It 1s controlled by subroutine SIMGO 
[Ref. 2). 


The Fortran77 code for the simulation is as follows: 


PROGRAM Mass Attack with Invisible Kill 


REAL ProbAcq, TimeWindow, TimeAcq, TimetoKill, Rate, Average, 
- Variance, Interval Width, Lower, Upper 

INTEGER Iseed, NrMissiles, m, Processed, NrAcquired, k, 
+MissilesDown, 1 


COMMON/ SHARE2/ Processed, NrAcquired, MissilesDown 
COMMON/ SHARE3/ Iseed, ProbAcq, TimeAcq, TimetoKill, Rate, 
+NrMissiles 

COMMON/ SHARE4/ Time Window 

COMMON/ SHARES/ Average, Variance, IntervalWidth, Lower, Upper 


Create input and output files. 
OPEN( 15, FILEZ'/EXPINV DATA) 
OPEN( 20, FILEZ'/EXPINV OUTPUT”) 


Initialize the seed for the random number generator 
READ(15, *) Iseed 
По 1 = 1, 16 
Read one combination of input parameters (see Table 5.1), and replicate the 
simulation for each of them. 
Read (15,*) Rate, NrMissiles, ProbAcq, TimeAcq, TimeWindow, m 
The time assigned to shoot at each missile is computed from the optimal 
number of missiles killed according to equation (2.2). 
If( m .ne. O ) then 
TimetoKill = TimeWindow / m - TimeAcq 
Else 
TimetoKill = 0.0 


46 


End If 


Initialize the statistics. 
Variance = 0.0 
Average = 0.0 
Interval Width = 1.0 


Perform independent replications of the simulation. 
Do k = 1, 15000 
CALL INIT 
CALL SIMGO 
CALL STAT( MissilesDown, k) 
End Do 
Write(20,*) 1, Lower, Average, Upper 
End Do 


CLOSE(15) 
CLOSE(20) 


Stop 

End 

SUBROUTINE INIT 
Initialize variables used in the simulation before each replication. 
REAL Now, Time Window 

INTEGER Processed, NrAcquired, MissilesDown 
COMMON/ SHARE1/ Now 

COMMON/ SHARE?2/ Processed, NrAcquired, MissilesDown 
COMMON/ SHARE4/ TimeWindow 

Now = 0.0 

Processed = 0 

NrAcquired = 0 

MissilesDown = O 

CALL CALANDAR INITIALIZATION 


Schedule the first acquisition event 
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CALL SCHEDULE( Now, 1) 


Schedule end of attack after time window elapsed 
CALL SCHEDULE( Now+Time Window, 3) 


Return 
End 
SUBROUTINE ACQUISITION 


REAL X, ProbAcq, Now, TimeAcq, TimetoKill, Rate 
INTEGER Iseed, NrAcquired, Processed, MissilesDown, NrMissiles 
LOGICAL Acquired 


COMMON/ SHAREI/ Now 

COMMON/ SHARE2/ Processed, NrAcquired, MissilesDown 
COMMON/ SHARE3/ Iseed, ProbAcq, TimeAcgq, TimetoKill, Rate, 
+NrMissiles 


SAVE SHARE2 


Sample from a Uniform(0,1), to determine if acquisition was successful 
CALL RANNUNM( 1, Iseed, 0.0, 1.0, 0.0, X) 
Acquired = .FALSE. 
If( X .le. ProbAcq ) Acquired = . TRUE. 
If( Acquired ) then 
NrAcquired = NrAcquired + 1 
Schedule a shooting event for the acquired missile, after acquisition time 
elapsed. 
CALL SCHEDULE( Now-TimeAcq, 2) 
Return 
else 
Processed = Processed + 1 
If( Processed .It. NrMissiles ) then 
If not acquired and not all missiles have been processed, schedule another 
acquisition. 
CALL SCHEDULE( Now+TimeAcqg+TimetoKill, 1) 
else 
End the simulation if all missiles have been processed 
CALL SCHEDULE( Now+0.00000000001, 3) 
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End If 
End If 


Return 
End 
SUBROUTINE SHOOTING 


REAL Now, TimeAcq, ProbAcq, Rate, Tsample, TimetoKill 
INTEGER MissilesDown, Processed, Iseed, NrAcquired, NrMissiles 
LOGICAL Destroyed 


COMMON/ SHAREI/ Now 

COMMON/ SHAREZ/ Processed, NrAcquired, MissilesDown 
COMMON/ SHARE3/ Iseed, ProbAcq, TimeAcq, TimetoKill, Rate, 
+NrMissiles 


SAVE SHARE2 


Samples the task time from an Exponential distribution 
CALL RANNUM( 3, Iseed, Rate, 0.0, 0.0, Tsample) 
Destroyed = .FALSE. 


If the random task time is less than the time assigned to shoot the missile, it 
will be destroyed. 
If( Tsample .le. TimetoKill ) Destroyed = .TRUE. 
If( Destroyed ) MissilesDown = MissilesDown + 1 
Processed = Processed + 1 
If( Processed .lt. NrMissiles ) then 

Schedule another acquisition event 

CALL SCHEDULE( Now+TimetoKill, 1) 
else 

if all missiles processed, end the attack 

CALL SCHEDULE( Now, 3) 
End If 


Return 


End 
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В. THRESHOLD MODEL 

Most of what was written for the Invisible Kill case still applies. All variables 
retain their meanings. The variable "Threshold" is introduced to denote the optimal 
threshold, «x^, which must be supplied to the simulation. In the current case 
optimal thresholds were computed by numerical solution of (3.10). 

Fortran listing follows: 

PROGRAM Threshold Simulation 


REAL ProbAcq, TimeWindow, TimeAcq, Threshold, Rate 
t, Average, Interval Width, Variance, Lower, Upper 
INTEGER Iseed, NrMissiles, Replications, MissilesDown, k 


COMMON/ SHARE2/ MissilesDown 

COMMON/ SHAREJ/ Iseed, ProbAcq, TimeAcq, Threshold, Rate, 
+NrMissiles 

COMMON/ SHARE4/ Time Window 

COMMON/ SHARES/ Average, Variance, IntervalWidth, Upper, Lower 


OPEN( 15, FILE=’/THRESH DATA’) 
OPEN( 20, FILE=’/THRESH OUTPUT’) 


Read(15,*) Iseed 
Do k= 1, © 
E Read input parameters and optimal threshold 
Read(15,*) Rate, NrMissiles, ProbAcq, TimeAcq, TimeWindow, 
+ Threshold 


Variance = 0.0 
Average = 0.0 
Replications = 0 
Interval Width = 1.0 


$ Replicates at least 10,000 times and until 95% confidence interval is less than 
т or equal to 5% of the average. 
Do While( IntervalWidth .gt. 0.05*Average .or. 
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+ Replications .le. 10000 ) 
Replications = Replications + | 
CALL INIT 
CALL SIMGO 
CALL STAT( MissilesDown, Replications) 
End Do 
Write(20,*) k, Lower, Average, Upper 
End Do 


CLOSE(15) 
CLOSE(20) 


Stop 
End 
SUBROUTINE ACQUISITION 


REAL X, Now, ProbAcq, TimeAcq, Threshold, Rate 
INTEGER Iseed, MissilesDown, NrMissiles 
LOGICAL Acquired 


COMMON/ SHAREI/ Now 

COMMON/ SHARE?2/ MissilesDown 

COMMON/ SHAREJ23/ Iseed, ProbAcq, TimeAcq, Threshold, Rate, 
+NrMissiles 


Sample from a Uniform(0.1) to determine outcome of acquisition. 
CALL RANNUNM( 1, Iseed, 0.0, 1.0, 0, X) 
Acquired = .FALSE. 
If( X .le. ProbAcq ) Acquired = . TRUE. 
If( Acquired ) then 
When missile is acquired schedule a "Shooting" event 
CALL SCHEDULE( Now-TimeAcq, 2) 
else 
Upon failure to acquire schedule another acquisition in time = current time + 
acquisition time + threshold. 
CALL SCHEDULE( Now+TimeAcq+Threshold, 1) 
End If 


Return 


s 


End 
SUBROUTINE SHOOTING 


REAL Now, Tsample, Threshold, Rate, ProbAcq, TimeAcq 
INTEGER MissilesDown, Iseed, NrMissiles 
LOGICAL Destroyed 


COMMON/ SHARE1/ Now 

COMMON/ SHARE2/ MissilesDown 

COMMON/ SHARE3/ Iseed, ProbAcq, TimeAcq, Threshold, Rate, 
+NrMissiles 


SAVE SHARE2 


Sample the exponential destruction time 
CALL RANNUM( 3, Iseed, Rate, 0.0, 0.0, Tsample) 
Destroyed = .FALSE. 
When the sampled time is less than the threshold time, the missile is killed 
If( Tsample .le. Threshold ) Destroyed = .TRUE. 
If( Destroyed ) then 
MissilesDown = MissilesDown + 1 
If all missile have been shot ... 
If(MissilesDown .eq. NrMissiles) then 
Return 
else 
If some missiles remain, schedule another acquisition immediately 
CALL SCHEDULE( Now+Tsample, 1) 
End If 
else 
If the missile is not destroyed, schedule an acquisition in time = current time 
+ threshold time 
CALL SCHEDULE( Now+Threshold, 1) 
End If 


Return 


End 
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С. PERFECT TASK INFORMATION 
We present programs for the simulation, Perfect Task Info Simulation”, and 
for the probability model, "Exact PMF for Perfect Task Info”. This simulation, 
unlike the preceding, has no discrete event subroutines. It represents the case 
without the threshold improvement mentioned in Chapter V. The improved case, 
however, is simulated with minor modifications to the code presented here. The 
same apples to the plots displayed in Chapter IV, for the PMF of N(T). To generate 
those plots, we also used a small APL function which computes the frequency of 
the data, adapted from Thomson [Ref. 10]. 
The main program variables are: 
e "E" is a real array for storage of the task times. 
e "Sum" accumulates the expended time. 
e “Time” is the time window, T. 
* "Rate" is the rate parameter, A, of the exponential task times. 


„и 


e "р" 15 the probability of success in acquisition. 

* "w' is the fixed setup or acquisition time. 

e "S" is used to store samples from a Uniform distribution. 
e "k", "kk", "kkk", are used as indices for "Do" loops. 


м 


eo "n" is the number of missiles. 
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The commented simulation program is listed next. Subroutine SORT 
implements an insertion sorting algorithm, which closely follows the code presented 
by Etter (Ref. 11]. Subroutine STAT is also used by the previous simulations. 

PROGRAM Perfect Task Info Simulation 
Simulation for the improved case 


REAL E(500), Sum, Time, Rate, p, w, S, Lower, Average, Upper, 
+ Variance 
INTEGER Iseed, Replications, k, kk, kkk, n 


COMMON/ MAINSTAT/ Lower, Average, Upper, Variance 


OPEN( 10, FILE=’/PACKGT DATA’) 
OPEN( 15, FILE=’/PACKGT OUTPUT’) 


Read(10,*) Iseed, Replications 


Do kkk = 1, 16 
Read(10,*) Rate, n, p, w, Time 
Average= 0.0 
Variance = 0.0 


Do kk = 1, Replications 
* Sample the n task times from an Exponential 
Dokzl,n 
CALL RANNUM( 3, Iseed, Rate, 0.0, 0.0, E(k) ) 
End Do 
$ Sorts the task times in ascending order of their lengths 
CALL SORT( n, E) 
k=0 
Sum = 0.0 
Do While( k .le. n .and. Sum .lt. Time ) 
k=k+1 
x Sample from a Uniform to determine if acquisition successful 
CALL RANNUM( 1, Iseed, 0.0, 1.0, 0.0, S) 
ЦС .le. p) then 
* If acquired pack one more job 
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X X * ++ 


Sum = Sum + w + E(k) 
Else 
If not acquired, expend time window 
Sum = Time + 1.0 
End If 
End Do 
Update the statistics 
CALL STAT(k-1,kk) 
End Do 


Write(15,*) kkk, Lower, Average, Upper 
End Do 


CLOSE(10) 
CLOSE(15) 


Stop 
End 
SUBROUTINE SORT( n, E) 


Sorts the elements of the array E in increasing order. E(1) = smallest; 
E(n) = largest. 


INTEGER n, kk, k 
REAL E(n), Stor 
LOGICAL Done 


Do kk = 1, n-1 
If( E(kk) .gt. E(kk+1) ) then 
Done = .FALSE. 
k = kk 
Do While( .not. Done ) 
Stor = E(k) 
E(k) = E(k+1) 
E(k+1) = Stor 
If( k .eq. 1 .or. E(k) .gt. E(k-1)) then 
Done = . TRUE. 
else 
k = k-l 


29 


End If 
End Do 
End If 
End Do 


Return 

End 

SUBROUTINE STAT({( y, n) 
Computes sequential statistics 


INTEGER y, n 
REAL Average, Variance, Oldaverage, Upper, Interval width, 
+Oldvariance, Lower, X 


COMMON/ MAINSTAT/ Lower, Average, Upper, Variance 
SAVE MAINSTAT 


Oldaverage = Average 
Oldvariance = Variance 
X = Real(y) 


If( n .eq. 1 ) then 
Average = X 
Variance = 0.0 
Intervalwidth = 1.0E20 
End If 


If(n .ge. 2) Average = (REAL(n-1)*Oldaverage+X)/Real(n) 
If(n .eq. 2) Variance - Oldaverage* *24(X**2)-(2.0*(Average**2)) 
If(n .ge. 3) then 
Variance = (REAL(n-2)/REAL(n-1)) * Oldvariance + 
+ Oldaverage**2 - (REAL(n)/REAL(n-1)) * Average**2 + 
+ (X**2)/REAL(n-1) 
End If 
Іп .GE. 2) Intervalwidth = 2.0*1.96*SQRT(Variance)/SORT(REAL(n)) 
Lower = Average - Intervalwidth/2.0 
Upper = Average + Intervalwidth/2.0 
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Кецип 


End 


The program “Exact PMF for Perfect Task Info" implements the analytical 


results derived in Chapter IV. It should be run with option "(NOXUFLOW’, to 


avoid underflow messages. 
Additional variables in this program are: 
e "PMF" - represents P(N(T) =). 
* "FSj" - denotes P(S,>T). 
The remaining variables retain their meanings and notation. 
PROGRAM Exact PMF for Perfect Task Info 


DOUBLE PRECISION Rate, Time, p, w, PMF, FSj, Average 
INTEGER j, n, k 


Ореп( 20, File 2 '/PEX DATA") 
Open( 25, File = '/PEX OUTPUT?) 


Do k = 1, 16 
Read(20,*) Rate, n, p, w, Time 
Average = 0.0D0 
Оо } = 0, п 
P Uses equation (4.23) to compute PMF 
PMF = (p**j)*(1.0D00-FSj(j,Time-DBLE(])*w,n,Rate)) - 
+ (p**(j+1))*(1.0D00-FSjG+1,Time-DBLE(Gj+1)*w,n,Rate)) 
E Average is computed with (5.1) 
Average = Average + DBLE() * PMF 
End Do 
Write(25,*) k, Average 
End Do 
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Close(20) 
Close(25) 


Stop 
End 
DOUBLE PRECISION FUNCTION FSjGj,X,nn,Rrate) 


Implements equation (4.17) 
INTEGER jj, i, k 
DOUBLE PRECISION X, Sum, Prod, Rrate, COMB, FACT 


If( X Je. O) then 
FS; = 1.0D00 
Return 
End If 
If( jj le. 0) then 
FS} = 0.0D00 
Elself((j).ge.1) .and. (j).le.(nn-1))) then 
Sum = 0.0000 
Compute the product term in equation (4.13) 
Doi = |, jj 
Prod = 1.0D00 
If ( jJ .ne. 1 ) then 
Do k = 1,3) 
If ( 1 .ne. k ) then 
Prod=Prod/(DBLE(nn-k+1)/DBLEQGj-k+1)-DBLE(nn-i+1)/ 
+ DBLE(j-1+1)) 
End If 
End Do 
End If 
I(ABS(Rrate*(DBLE(nn-i+1)/DBLE(j-i+1))*X) .lt. 170.0) then 
Sum=Sum+Prod*(DBLE(jj-i+1)/DBLE(nn-i+1))*DEXP(-Rrate* 
(DBLE(nn-i+1)/DBLE(jj-i+1))*X) 
End If 
FSj = Sum * COMB(nn,jj) 
If( ABS(FSj) .gt. 1.0D00) FS; 2 1.0000 
End Do 
Elself(jj .eq. nn) then 
Sum = 1.0D00 
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Do K = nn) 
Sum = Sum + ((Rrate*X)**k)/FACT(k) 
End Do 
FS; 2 Sum * DEXP(-Rrate* X) 
Elself(jj .gt. nn) then 
FS) = 1.0D00 
End If 


Return 
End 
DOUBLE PRECISION FUNCTION COMB( m, p) 
Calculates the binomial coefficients ы 
р 

INTEGER m, p, k 
COMB = 1.0D+00 
If( p .ne. m) then 

Dok=1,p 

COMB = COMB*(DBLE(m-k+1)/DBLE(k)) 

End Do 
End If 
Return 
End 


DOUBLE PRECISION FUNCTION FACT(X) 


Returns the value of X! 
INTEGER X, K 


FACT = 1.0D+00 
If( X .ge. 2 ) then 
"Оо К = 2, Х 
FACT = FACT * DBLE(K) 
End Do 
End If 
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Return 


End 
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APPENDIX B. THRESHOLD MODEL GRAPHICS 


This appendix contains plots for the Threshold model. The response surfaces 
shown were obtained in GRAFSTAT [Ref. 12], from simulation output. 
In both plots the most salient feature is that the expected number of missiles 


killed is not very sensitive to the threshold alone. 
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Figure Bl. Threshold and Acquisition Time as factors 
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E[ missiles killed] 


Figure B2. 
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APPENDIX C. DISTRIBUTION OF S, 


In this appendix, we describe the evidence which suggested that a Gamma 
distribution would provide a suitable approximation to the distribution of S, (the 
approximation is developed in subsection IV.C.1). 

A slight adaptation of the program to simulate the packing problem (Perfect 
Task Information) was used, to produce samples of S; for several values of j. These 
samples of size 1000, were used to construct Gamma probability plots. One such 
plot is depicted in Figure Cl, for j=5. 

Other plots, with different values of j, produced much similar results. As can 
be seen, the sampled distribution seems to closely follow the shape of a Gamma. In 


all plots slight departures from linearity were only observed at the right tail. 
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GAMMA PROBABILITY PLOT, N=1000 
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